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ABSTRACT 

We present an investigation of sample selection effects that influence the observed black hole - bulge relations and its evolution with 
redshift. We provide a common framework in which all kinds of selection effects on the M, -bulge relations can be investigated, but 
our main emphasis is on the consequences of using broad-line AGN with virial estimates of black hole masses and their host galaxies 
to search for evolution in the BH - bulge relation. We identified relevant sources of bias that were not discussed in the literature so 
far. A particularly important effect is caused by the fact that the active fraction among SMBHs varies considerably with BH mass, 
in the sense that high-mass BHs are less likely to be active than lower mass ones. In the connection with intrinsic scatter of the 
BH - bulge relation this effect implies a bias towards a low BH mass at given bulge property. This effect adds to the bias caused 
by working with luminosity or flux limited samples that were already discussed by others. A quantitative prediction of these biases 
requires (i) a realistic model of the sample selection function, and (ii) knowledge of relevant underlying distribution functions: The 
distribution of spheroid properties such as velocity dispersions or masses; the active fraction as a function of BH mass, or alternatively 
the active black hole mass function; and the Eddington ratio distribution function. We employed our formalism together with recently 
determined distribution functions to investigate how much existing studies of the M.-bulge relation using AGN hosts might suffer 
from selection biases. For low-redshift AGN samples we can naturally reproduce the flattening of the relation observed in some 
studies without having to invoke intrinsic differences in the BH - bulge relation between active and inactive galaxies. When extending 
our analysis to higher redshift samples we are clearly hampered by limited empirical constraints on the various relevant distribution 
functions. Using a best-guess approach for these distributions we estimate the expected magnitude of sample selection biases for 
a number of recent observational attempts to study the M.-bulge evolution. In no case do we find statistically significant evidence 
for an evolving M.-bulge relation. While the observed apparent offsets in M. /MB U i ge from the local relation can be quite large, the 
sample selection bias estimated from our formalism is typically of the same magnitude. We suggest a possible practical approach to 
circumvent several of the most problematic issues connected with AGN selection; this could become a powerful diagnostic in future 
investigations. 
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1. Introduction 

Supermassive black holes a ppear to be ubiquitous in th e 
centers of massive galaxies dKormendv & RichstoneL I1995I) . 
Furthermore, they show tight correlations with the proper- 
ties of the spheroidal components of their host galaxies, e.g . 
with the stellar velocit y dispersions ( Ferrarese & Merrittl 200C ; 



not well established dGiiltekin et all I2009I) . Furthermore, there 
is evidence that late-type galaxies and in particular pseudob- 
ulges do not follow the same relations as early-type galaxies 



I2011 



2008r. iGrahaml 2008} iGreene et all l2010t iKormendv etafl 



Gebhardt et all l2000t; iTremaine et all l2002t iGultekin et al 



1998 



20091). bulge luminosit i es and masses jMagorrian etaL 
Marconi & Hunti 120031; fearing & Rixl 12004 ISani et all 1201 ll) 
or concentration indices dGraham et all l200lb . These rela- 
tions have been established by direct dynamical measurement 
of the black hole mass in a few doze ns of nearby galax- 
ies, mainly using stel l ar dyn a mics (e.g. [van der Marel et all 
Il998t lEmsellem et all Il999t [Gebhardt et all l2003h or~gas 



Jiang et alll201ll) . 
The black hole - bulge relations contain important in- 
formation, as they imply a connection between the growth 
of a black hole and the evolution of its host galaxy. A 
common picture invokes AGN feedback to shut down star 
formation and self-regulate black hole accretion. Models 
building on this scenario are able to reproduce the lo- 
cal black hole - bulge relations in numerical simulations 



dynamics (e.g. | Ferrarese et al 1 119961 iMarconi et al 1 120011 
Ml 



iPalla Bonta et all 120091) 



dDi Matteo et all 120051: ISiiacki et all 120071: iBooth & Schave . 



20Q9Tr and semi-anal ytically dKauffmann & Haehneltl 



2000; 



2006 



Nevertheless, significant uncertainties in these relations per- 
sist. The high mass and the low mass regimes are still poorly 
constrained. There may also be systematic effects on the de- 
termined masses; for example, not accountin g for triaxiality of 
the galaxy (Ivan den Bosch & de Zeeuwll2010l) or not including a 
dark matter halo in the dynamical models dGebhardt & Thomas! 
120091: ISchuTze & GebhardA I20TTT) can bias the dynamical mea- 
surement. Also the intrinsic scatter in these relations is still 



Cattaneo et all 120051: [Croton et all 120061 [ Bower et all 
ISomerville et all 120081: iMarulli et all 120081) . A correlation be- 
tween black hole mass and host galaxy mass will also be tight- 
ened or may even be generated natur a lly within a merger driven 
galaxy evolution framework dPengl 120071 : iHirschmann et all 
l20lol:lJahnke & Macci6Ll20TTl) . 

Essential constraints on the origin of the black hole - 
bulge relation can be inferred from their redshift evolution. 
There are sever al theoretical pr e dictio n s, based on numeri- 
cal simulations dRobertson et all l2006t iHopkins et all 2007 
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Pi Matteo etail 120081: iJohansson et all l2009t IBooth & Schave , 



201 lb Dubois et all 1201 11) as well as on semi-analytic mod- 
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els (IWvithe & Loebl l2003t ICrotonL 120061: iMalbon et all l2Q07t 
iHopkins et all l2009tTLamastra et all l2010l) . Although the de- 



tails are still far from being settled, they tend to predict an in- 
crease in the M./Msaige ra tio with redshift, while only a weak or 
even negative evolution in the M.-cr, relation is expected. 

These models need to be confronted with observations. 
Several different approaches have been followed in the last 
years to observationally constrain the evolution in the black hole 
- bulge relations. These range from more indirect arguments, 
to direct estimates of M. and the respective bulge property. 
Constraints on integr ated quantities can be gained from the black 
hole mass function (iMerloni et all l2004t IHopkins et all p006: 
Shankar et all 120091: Somervilld 120091: iKisaka & Koiimal boiol 
Willottet all £01 Ob . Bluck et all (1201 lb studied X-rav selected 



AGN, employing Eddington ratio arguments, to constrain black 
hole - bulge coevolution. 

Direct dynamical determinations of M, are not feasible out- 
side of the very local universe. Therefore, the most commonly 
adopted approach resorts to broad line AGN, for which black 
hole masses ar e easily accessible through the so-called 'virial es- 
timato r' (e.g., iMcLure & Jarvisl 120021: IVestergaard & Petersonl 
12003) . The main challenge for the determination of the M. -bulge 
relationship from AGN samples is the determination of the bulge 
properties, hampered by the bright nuclear point source of the 
AGN, which may well outshine the entire galaxy. 

Evolution in the M.-<x„ relation has been studied either by 
meas u ring stellar velocity dispersions directly (IWoo et alll2006l 
l2008t IShen et all l20Q8h or by usi ng the widths of narrow 
emission lines as surro gates for cr, dShields et all 120031 l2006t 
ISalviander et all [20071) . The M.-M^ a \ s& rel ation has been inve s- 
tigated fro m QSO h ost galaxy luminositie s dPeng et all 2006allb ; 



Treu et all 120071: iMcLeod & Bechtoldl 12009b iDecarlietal 

et_all 

the stellar masses utilisin g colour information (ISchramm et 



2010t iBennert et all 1201 Ol iTargett et ail 1201 ll). by estimatin 



selection effects on observations of the M.-MB U i ge relation and 
its evolution. We also investigate how well these biases can be 
quantified and corrected for. The paper is organised as follows: 
In section [2] we present our general framework to investigate 
the consequences of selection effects and show a first applica- 
tion to the quiescent black hole sample. Section[3]discusses sev- 
eral selection effects for AGN samples and their ramifications 
for the black hole-bulge relations. In section [4] we take into ac- 
count redshift evolution effects. We discuss the implications of 
our results on observational studies in section[5] We finally con- 
clude in section [6] For the cosmological parameters we assume 
H = 70 km s- 1 Mpc" 1 , fi m = 0.3 and Q A = 0.7. 



2. The local M.-M Bu i g e relation 
2. 1 . General framework 

We start with some basic thoughts on the BH - bulge relation 
and its sensitivity to selection effects, postponing specific aspects 
of AGN samples to the next section. We adopt the following 
convention: Black hole masses are given by fi = log M., while 
any of the various relevant bulge properties is denoted as s, with 
s = log cr, or s = log MBuiee, respectively. In some examples we 
focus on the M.-Ms u i ge relation, but the arguments equally hold 
for the M.-cr, relation. We explicitly do not imply one relation 
to be more fundamental than the other. 

The distribution of objects in the jj.-s diagram is given by the 
bivariate distribution function of bulge property and black hole 
masses VF(s,/i). Thus y ¥(s,/j)dfids gives the number of objects 
per Mpc 3 with galaxy property between s and s + ds and black 
hole mass between jj and y. + d/i. We assume the following pa- 
rameterisation for most of the paper, 



ing 

2 y(s,fi)=g(n\s)® s (s), 



(1) 



1 



V2^c 



■ exp < - 



I2008t Uahnke et all 120091: Merloni et all [201 Ob . or by dynamical 
mass measurements (Ilnskipetal.r2011 ). A lso obscured AGN 
with detectable broad lines have been used (ISarria et al 
iNesvadba et all 1201 lb . for which the determinat ion of stellar 
masses is less problematic. IMcLure et al.l (12006b inferred the 
M.-Mfluige ratio up to z — 2 by matching distributions of radio- 
loud QSOs and radio galaxies. At the highest redshifts, CO rota- 
tion curves have been used to determine dynamical masses for a 
few individual obj ects dWalter et all 12004b iRiechers et all 120081 
120091: IWang et alll2010b . All these methods have their own ad- 
vantages and drawbacks. Nevertheless, there is increasing evi- 
dence for apparent evolution, in a sense of high-z SMBHs being 
more massive at given bulge mass than in the local universe. 

However, it has been suspected that these results may be 
biased due to sample selection effects. One such bias could 
arise from the fact that luminous AGN on average tend to 
have massive black holes. The intrinsic scatter in the BH - 
bulge relation then generates a bias towards a higher M./MsuUe 
ratio. Such a bias has been pointed out by several authors a / \ _ fm/ ~. H _ a / \ 
dAdelberger & Steidelll200l:lFine et alll2006tlSalviander et all s{s) ~ J T(i 'W d ^ ~ ® M 



where O s (s) is the spheroid distribution function, e.g. the 
spheroid mass function, and g(fi \ s) gives the probability of find- 
ing the black hole mass yu, given s. If s and fi are correlated, 
as suggested by the observations, then g(/i \ s) corresponds to 
this correlation. In the following we assume a linear relation 
H — a + bs with log-normal intrinsic scatter cr, i.e. 



(li - a - bs) 2 
~2a 2 



(2) 



This parameterisation for the bivariate distribution function is 
not the only one possible. We comment on the observational jus- 
tification of this form in the Appendix. While our quantitative re- 
sults depend on the adopted parameterisation, most of our qual- 
itative results are independent of the specific choice of ^(s,^. 

The corresponding distribution functions for s and are 
given by marginalisation over the other variable, 



but received major attention by the work of Lauer et all 
I who discussed it in some detail ( but see also |PengLl2010b . 
A further effect recently discussed by IShen & Kellvl d2010b is 
caused by the uncertainty in virial black hole mass estimates in 
connection with the steep decrease of the active black hole mass 
function. 

In this paper we argue that besides the above-mentioned ef- 
fects there are even more potential biases that need to be consid- 
ered. While previous authors mostly focused on one particular 
source of bias, we aim at providing a comprehensive study of 



O.0u) = j V(s,n)ds= j *(/i|j)<D,(j)dj 



(3) 



(4) 



The integration over fi simply returns the spheroid mass func- 
tion. The integration over s provides the quiescent black hole 
mass function (BHMF), equivalent to the common approach 
for its determin ation (e.g. lYu & Lul 120041: iMarconi et all 120041: 
iMerlonil l2004b . The first equality in Equations Q and © is in- 
dependent of the special choice of ^(i,//). In particular, as the 
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Fig. 1. Bivariate probability distribution function of bulge 
masses and black hole masses with no selection effects. The 
contours indicate decreasing probability in logarithmic units. 
The thick bla c k line shows the input M,-M^ u \„ e relation from 
lHaring & Rixl (120041) . The red dot-dashed line shows the mean 
galaxy property at a given black hole mass. The upper panel 
shows the projection of the bivariate distribution function to 
the bulge mass, i.e. the spheroid mass function. The right panel 
shows the projection to black hole mass, i.e the quiescent black 
hole mass function. 

galaxy distribution function is an observable, Equation (0 sets a 
constraint on the bivariate distribution y(s,/j). Note that the spe- 
cific shape of ^(s./i) has direct consequences for the quiescent 
black hole mass function. 

Fig.[T]shows *TO,yL(), 0. s (i) and <£>.(//) for the local universe, 
employing our parameterisation. We estimated the spheroid 
mass functi on from the earl y type and late type stellar mass 
functions of iBell et al.l (120031) . We computed the spheroid mass 
function as sum of the early type mass function and the mass 
function of the bulge components of late type galaxies. Since 
the bulge fraction in disk galaxies is still somewhat uncertain, 
we simply assumed an average value of B/T = 0.3 to convert 
the late type mass function into a b ulge mass function, con- 
sistent with cu rrent observations (e.g. lGraham & WorlevL l2008t 
iGadottil 12009). This value also produces a black hole mass func- 
tion that is con sistent with determin ations from the M.-cr, and 
M.-L relations dMarconi et al.ll2004j). For the M .-M Bu i ee relation 
we used the relation by lHaring & Rixl d2004l) and an intrinsic 
scatter of <x = 0.3 dex. For most of this paper we keep this value 
for the intrinsic scatter fixed; we comment on the consequences 
of changing it when appropriate. 

The distribution shown in Fig. [T] would be obtained for 
a strictly volume limited sample of galaxies with Meuige and 
M. measurements. In practice, the distribution will be affected 
by the way that the sample was constructed. This is formally 
accounted for by a selection function Q, defined as the proba- 
bility of observing an object of a given bulge mass, black hole 
mass and further possible selection criteria, such as for instance 
redshift or AGN luminosity. The selection function can thus be 



written as 0(s, where 6 refers to the set of additional pa- 
rameters present as selection criteria. The observed apparent bi- 
variate distribution is then given by 



-/ 



(5) 



where pg(ff) is a set of normalised distribution functions of the 
parameters 9. Depending on Q, the observed M.-Ms u i S e relation 
can be significantly biased if the selection effects are not taken 
into account. However, this requires a proper knowledge of the 
selection function. 

The bivariate distribution ^(s,//), when normalised to one, 
represents the full probability distribution of the expected M, - 
A^Buige relation. However, it may happen that we are not so much 
interested in the full distribution, but rather in the mean relation. 
This is given by the mean black hole mass at a given bulge prop- 
erty 



0i>(*) = 



J )u¥ 0,ju)dju 



(6) 



If no selection effects are present then (/i)(s) will be identical to 
the input M, -bulge relation, (/j) - a + bs. Alternatively we can 
consider the mean bulge property at a given black hole mass 



<*>(A0 



JV (s,fi)ds 



(7) 



Even without any selection effects, (s)Qi) does not simply corre- 
spond to the inverse relation. It deviates from it quite strongly 
at the high mass end where the space density in the galaxy 
distribution function is decreasing. This effect is illustrated as 
the red dot-dashed line in Fig. Q] and was already discussed by 
lLauer et all d2007l) . It is not an observational bias but just a direct 
consequence of the different projections chosen, combined with 
the steep decrease of the galaxy distribution function. In gen- 
eral, under the presence of selection effects {s)(jj) may change 
as well. 

In observational studies of the evolution of the M. -bulge re- 
lation, the sample is often compared to the local relation in terms 
of a single value, the offset from the local relation. The sample 
bias for this offset is given by 



ff(/d - a- bs) YqCs, /j.) dju ds 



(8) 



Alternatively, the mean offset of the galaxy property from the 
local relation is given by 



(As) = = --(Afi) , 

ffW (s,M)dfxds b 



(9) 



where we use the inverse relation s — a + B/u with a = -a lb 
and B = l/b. This directly leads to (As) = -l/b(A[i), without 
any specific assumption on ^(s,^). Thus the biases in (Afi) and 
(As) should always be directly proportional to each other. 

Already the well-studied samples of a few dozen local galax- 
ies with reasonably well-measured dynamically-based black 
hole masses will invariably be affected by selection effects. 
However, these samples are very inhomogeneous, and presum- 
ably they cannot be described by a single well-defined selec- 
tion function, which hampers a proper inv estig ation of poten- 
tial selection effects. lYu & Tremaind J2002) and iBernardi et al.l 
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Fig. 2. The bivariate probability distribution function and its projections for a sample of black holes with dynamically measured 
masses, under the presence of a sample selection depending on the resolution of the black hole's sphere of influence. A threshold 
of ^infMes > 1 (left panel) or R m f/d res > 0.1 (right panel) is applied to the sample. The dashed black line shows the input M.- 



A^Buige relation, the thick black solid line highlights the 
'observed' relation (s)(pi) 



observed' relation (jj)(s), and the red dashed dotted line indicates the 



d2007l) reported on a possible bias in the dynamical black 
hole mass sample identified through finding a discrepancy in 
the cr, - L relation compared to the SDSS. More recently, 
IGultekin et all d2009b discussed a bias induced by culling the 
sample based on resolving the spheres of influence of the black 
holes. Before discussing selection effects on AGN samples in de- 
tail, we investigate this bias to illustrate the general applicability 
of Equation (|5j- 

2.2. Bias by sphere of influence resolution 

Dynamical black hole mass determinations need to spatially re- 
solve the region of gravitational influence of the central black 
hole on the stellar velocity distribution. The size of this region 
is commonly estimated by the 'sphere of influence' of radius 
Rmf = GM.cr~ 2 . It has been argued that dynamical black hole 
masses are unreliable for Rm/cI K s < 1 and therefore sho uld be 
excluded from the sample (e.g. iFerrarese" & Fordl l2005t). Here 
af res is the spatial resolution of the observations. IGiiltekin et al.l 
d2009l) demonstrated that rejecting galaxies below a line M. oc 
cr t leads to a bias in the M.-bulge relations. iBatcheldorl d2010l) 
pointed out that it would even be possible to artifi cially generate 
aM.-c r, correlation via this selection effect, while lGiiltekin et all 
d2011l) excluded thi s possi bility at least for the most extreme 
case. IGiiltekin et all (120091) argued against culling the sample 
based on the sphere of influence to avoid this bias. However, 
even without any active rejection the same effect may be ef- 
fectively present, either due to an implicit target selection or 
through the mere ability to detect a black hole. The condition 
Rinf/dies > 1 is no strict limit for the reliability of a black hole 
detection, but it is clear that for Ri n f/d KS < 1 no black hole can 
be detected dynamically. 

We now illustrate the results of this bias. We assume that we 
can detect all black holes above a certain threshold in resolution 
of the sphere of influence, /?i n f/c/ r es > '"min- We keep the angular 



resolution of our hypothetical survey fixed at 0.1". The selection 
function is then given by 



. 1 for R m fld KS > r n 
Ws,*d)= o else 



(10) 



where d is the distance of the galaxy. The observed bivariate 
probability distribution follows as 



r 



Q.(s,fi,d) x ¥(s,fi)d 2 dd. 



(ID 



For the purpose of illustration we again use the M.- 
MB U ige diagram. We convert bulge masses into velocity disper- 
sions assuming a fixed mass-to -light ratio for the r-band and 
using the cr t - L r relation from iBernardil d2007h . We cover the 
distance range from 1 to 30 Mpc. In Fig. [2] we show the M.- 
A^Buige diagram for two thresholds in the sphere of influence res- 
olution, r m ; n = 1 and r m ; n =0.1. We additionally plot the 'ob- 
served' mean M.-Ms u \ ge relation {/-i)(s) for both cases as thick 
black solid line and the 'observed' mean relation (s){/j) as red 
dashed-dotted line. The bias induced by whether or not resolv- 
ing the sphere of influence leads to a flatter slope, a higher 
no rmalisation and a sma ller intrinsic scatter, as already shown 
by IGultekin et al.l d2009l) . See also their work for an extensive 
Monte Carlo investigation of this bias. 

While a significant bias is induced only for a high threshold 
r m j n , the probability distribution is affected even for lower val- 
ues, most severely at the low mass end. Also the relation (s)(jj) 
is biased in both cases. The selection along lines of M. oc cr\ 
is clearly visible in the probability contours. However, it is hard 
to quantify how much the actually observed samples of black 
holes are actually biased by this selection effect, due to their 
inhomogeneous character. But these calculations show that se- 
lection effects are a concern in general, also for dynamical black 
hole mass measurements. 
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3. Biases of broad-line AGN samples 

Selection effects are a major concern for samples selected as 
broad line AGN. For these objects the black hole mass can be 
easily estimated using the virial method. This makes them so 
valuable for studies of the evolution in the M. -bulge relations, 
as at high redshift black hole masses cannot be determined by 
means of direct dynamical observations. We first discuss generic 
issues that are also inherent in low-z AGN samples. The appli- 
cation to higher redshifts will be covered in section|4] 

Samples of non-active galaxies are typically drawn from the 
galaxy luminosity function. In contrast, AGN samples are effec- 
tively drawn from the luminosity function of broad-line AGN. 
This can induce non-t rivial selection effe cts on the sample, as 
already emphasized bv lLauer et al.l (120071) . Firstly, there will be 
a luminosity bias. Flux limited AGN surveys will on average 
contain more luminous objects then an 'ideal' volume-limited 
sample. Thus, the AGN will on average host more massive black 
holes than in the volume-limited case. The intrinsic scatter in the 
black hole-bulge relations, together with the steep decrease of 
the spheroid distribution function, will then cause a Malmquist 
type bias towards more massive black holes at a given spheroid 
mass. The strength of the bias depends on the luminosity limit, 
the shape of the spheroid distribution function, the intrinsic scat- 
ter in the M, -bulge relation, and the distribution of Eddington 
ratios. We will further discuss and illustrate these details below. 

A second effect, so far not discussed in the literature, can 
be called an 'active fraction bias'. For any sample of broad line 
AGN, their black holes are not randomly drawn from the entire 
black hole population (as described by the BH mass function), 
because they are by definition active. Only a minority of all black 
holes is in an active state, quantified by the active fraction (or 
equivalently by the AGN duty cycle). If this active fraction is a 
function of M. then the intrinsic scatter in the M.-MB U i ge relation 
will induce an additional bias, as will be shown below. If, for 
example, the active fraction decreases with increasing black hole 
mass, then for a given spheroid mass it will be more probable to 
find a smaller mass black hole in an AGN sample, causing a bias 
towards a lower M./MB U i ge ratio. Conversely, for an increasing 
active fraction we expect a positive bias, while for a constant 
(M. -independent) active fraction no such bias will occur. Thus, 
this bias can work in both directions, depending on the black 
hole mass dependence of the active fraction. It will add to the 
luminosity bias. 

3. 1 . Luminosity limited samples 

We now show analytically how these selection effects influence 
the M. -bulge diagram. Broad line AGN are drawn from the AGN 
luminosity function <1>l(0, with I = logLboi- The luminosity of 
an AGN is produced by mass accretion onto a SMBH, and it is 
therefore determined by the black hole mass and the mass ac- 
cretion rate. The latter is expressed in normalised units by the 
Eddington ratio A = Lboi/^Edd- The Eddington luminosity is 
proportional to M. , thus the bolometric luminosity is given by 
I = log A + p + 38.1. The AGN luminosity function <& L (l) (in 
logarithmic units) is then given by 



-I 



p A (l - p) <D., a (//) dp . 



(12) 



By <S>. A (p.) we define the active BHMF, where active black holes 
are defined as those that contribute to the corresponding AGN 
luminosity function. Thus, when restricting the AGN sample to 
type 1 AGN (showing broad lines), the active BHMF includes 



only these type 1 AGN. The distribution function p^(l-p) gives 
the probability of finding a black hole with mass p given an 
AGN luminosity /, i.e. it corresponds to the normalised distri- 
bution function of Eddington ratios A. We implicitly assume 
that the Eddington ratio distribution function is independent 
of black hole mass, i.e. px(p \ I) = P a{1 - (i) = /^(log A). 
Type 1 AGN are observed to have Eddington ratios in the range 
0.01 < Lboi/^Edd < 1- When applying a lower threshold to the 
Eddington ratio, not all black holes are currently in an active 
state. The active fraction is defined as the ratio between the ac- 
tive and the total black hole population, estimated as 



p-Ap) = O., a 0u)/<5., q (ju) ; 



(13) 



with <t>. ia (yu) and <t>. q (ji/) as the active and quiescent BHMFs, 
respectively. Therefore, the AGN luminosity function can be ex- 
pressed as 



-JOT 



Ol(0 = PM - P) P-Ap) g(p I s) 0. s (s) dsdp . 



(14) 



There is thus a connection between the bulge property and the 
AGN luminosity, which is however smeared out by convolution 
with a set of additional distribution functions. 

Following Equation (0, the bivariate black hole-bulge distri- 
bution function for a luminosity limited sample is given by 

V (s,p) = Jn(s,pJ) PA (l-p)g(p\s)® s (s)dl. (15) 

If we assume a fixed lower luminosity limit /ij m , the selection 
function is 



Q.(s,p,l) 



dp) 




for / > Z lim 
else 



The bivariate distribution function is then 



*F (s,//) = p 3C (p)g(p\s) 



p x {l -p)dl. 



(16) 



(17) 



Hence, O/X0 and ^V {s,p) are just different projections of a mul- 
tivariate distribution function ¥ (s,p,l). The non-trivial differ- 
ence between ^(s^) and ^(s,^) is the source of the bias. It is 
controlled by the three probability distributions, p dc {p), p,\{l—p), 
and g(p | s). The distribution of pAp) regulates the active frac- 
tion bias, px(l - p) rules the luminosity bias and g(p \ s) controls 
the overall strength of the bias. 

In Equation ( fT7b we can combine the parts governing the 
AGN selection by defining the selection function, integrated over 
AGN luminosity, 



&(P) = PAp) Pa(1 ~P)dl. 



(18) 



and thus *F (s,//) = SlipfV^p). 

The mean relations obtained from such a sample, affected by 
an AGN luminosity limit, are 



ip){s) 



is)(p) 



J p V (i,//) dp j pQ.(p)g(p \s)dp 
fV (s,p)dp ~ j£l(p)g(p\ S )dp ' 

/ s x ¥ (s,p)ds _ j sg(p | s)<& s (s) ds 
fW (s,p)ds ~ fg(p\s)® s (s)ds ' 



(19) 



(20) 



The relation (p)(s) is independent of the galaxy distribution 
function, but it is affected by the selection function. Thus, 
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Fig. 3. Bivariate probability distribution function and its projections for a local type 1 AGN sample under the presence of a lower 
luminosity limit. The dashed black line again shows the input relation, the thick black solid line the observed relation (/u)(s), and the 
red dashed dotted line indicates the 'observed' relation (s){/j) of this sample. A lower luminosity limit of log L\, \ > 45 (left panel) 
and log Lboi > 47 (right panel) is applied, respectively. 



(l-i)(s) + a + bs, as would be the case without selection effects. 
On the other hand {s){/j) is already intrinsically affected by the 
galaxy distribution function, but it is independent of the selection 
function. Therefore, (s)(/i) contains information on the intrinsic 
relation, unaffected by the AGN luminosity bias and active frac- 
tion bias. This is an interesting fact, as it provides a potential 
route to study the evolution in the M. -bulge relations without 
having to account for AGN-specific selection effects. However, 
even if (s){/j) is not biased by the selection effects, the mean off- 
set (As) is still biased, as shown by Equation because we 
have to integrate over the entire observed black hole mass distri- 
bution. 

We now first illustrate the AGN bias in the M.-bulge plane, 
and then investigate it in terms of a simple offset from the input 
relation. For this purpose, we need to know all the underlying 
distribution functions, for which we here adopt the z ~ values. 
For $> s (s) we again use the above estimated spheroid mass func- 
tion. For g{fi | s) we assume a log-normal distributi on with intrin- 
sic sca tter cr = 0.3 dex around the relation from lHaring & Rixl 
d2004l) . These then define the quiescent BHMF via Equation ©. 
If we further know the active BHMF, the active fraction is given 
by Equation ( fT3l . The active BHMF and the Eddington ratio 
distribution function (ERDF) p,\{l — n) for type 1 AGN in the 
local universe were determined by ISchulze & Wisotzki (1201 Oi 
hereafter SW10), and we use those values, adopting a modified 
Schechter function for the BHMF and a Schechter function for 
the ERDF. The SW10 active BHMF implies a significant de- 
crease of the active fraction with increasing black hole mass, 
thus the active fraction bias will work towards a lower M, /Mfluige 
ratio, opposite to the luminosity bias. Note that this parameter- 
isation by construction ensures the consistency with the broad 
line AGN luminosity function via Equation (fT4l . 

In Fig. [3] we show the M,-M^ u \„ e diagram for a luminos- 
ity limited local AGN sample, for two lower luminosity limits, 
'min = 45 and l m i„ = 47. At z ~ the latter case is rather unrealis- 
tic, as QSOs of these luminosities are extremely rare. However, 
as further discussed below, it is of importance for higher z ob- 



servations, where such luminosities are quite common in QSO 
samples. A luminosity limit of l m \ n = 47 roughly corresponds to 
the SDSS magnitude limit of i ^ 19 at z - 3 and to the SDSS 
high-z limit of i = 20.2 at z - 6. 

For both luminosity limits the M.-Ms u \ge relation is biased 
towards larger black hole masses at a given spheroid mass. 
This is already true for Z m i n = 45, but it is greatly enhanced 
for the high luminosity limit. The solid line again shows the 
mean black hole mass for a given spheroid mass (computed via 
Equation ©), i.e. the apparent M.-Mn u i ge relation obtained for 
the given sample. It is clearly visible that in general the bias 
is stronger for lower masses. Specifically, it is strongest close 
to the luminosity limit of the survey, which on average corre- 
sponds to the lowest mass black holes. The red dashed-dotted 
line again shows the mean spheroid mass for a given black hole 
mass (computed via Equation (|7]i). As noted, this relation is un- 
affected by the AGN selection and is identical to the one shown 
in Fig. [1] However, it is not identical to the intrinsic relation, 
as it turns upwards at the high mass end due to the decrease in 
the spheroid mass function. Therefore, an integration over (s)(jj.) 
with a stronger weight on the high mass end will lead to a bias 
for the mean offset from the local relation, i.e for the total sample 
bias. 

The sample bias is greatly enhanced for a brighter luminosity 
limit. The steep decrease of the spheroid mass function and the 
intrinsic scatter in the relation itself induce an increase in the bias 
compared to the flat section of the spheroid mass function. This 
is shown in Fig. [4] where the mean offset from the input relation 
for the whole AGN sample, (A/i) (computed via Equation (JSJl), 
is plotted for a range of bolometric AGN luminosity limits, as- 
suming the local distribution functions. Here we assumed our 
reference value for the intrinsic scatter, cr = 0.3 dex, but we 
also indicate the dependence of the bias on the intrinsic scatter. 
A smaller intrinsic scatter also leads to a smaller sample bias, 
while a larger scatter causes a larger bia s. By the do t ted lin e we 
additionally show the bias predicted by lLauer et all (120071) . us- 
ing their Equation (25) and the local bolometric type 1 AGN 
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Fig. 4. Predicted overall bias for a luminosity limited local (z ss 
0) type 1 AGN sample as a function of the low luminosity limit 
(in logarithmic units). The solid line shows our prediction for 
cr = 0.3, while the dashed lines show the predictions for cr = 0.4 
(upper line) and cr - 0.2 (lower lin e). The dotted line shows the 
prediction from lLauer et alj (|2P07), but for our local AGN lumi- 
nosity function; see text for an explanation of the differences. 

luminosity function from lSchulze et al] d2009l) . The main differ- 
ence between their result and ours is that we not only use the 
AGN luminosity function, but evaluate the whole set of underly- 
ing distribution functions. In particular, we also take into account 
the active fraction bias that reduces the total bias, due to the de- 
crease of the active fraction towards higher black hole masses. 
As a caveat we note that for very high luminosities, our knowl- 
edge of the underlying distributions (active and passive BHMFs, 
AGN luminosity function and M.-M% u \ge relation) are observa- 
tionally poorly determined. 

We provide an additional illustration of these selection ef- 
fects by Monte Carlo Simulations. We generated a large galaxy 
sample drawn from the spheroid mass function, and attributed 
black hole ma sses drawn from a lo g-normal distribution with 
mean from the lHaring & Rixl (l2004h relation. We decided if the 
black hole is in an active state based on the probability p ac (ji). 
If it is active, an Eddington ratio was drawn from the SW10 
Eddington ratio distribution function, which also sets the bolo- 
metric luminosity of the AGN. Now a number of specific selec- 
tion criteria were applied to this sample, with the results pre- 
sented in Fig. |5] It is notable that for a bright luminosity limit 
most, if not all, AGN lie above the input M.-M% u \ se relation, 
while they are intrinsically drawn from it. The applied lumi- 
nosity limit effectively corresponds to an unsharp lower limit 
on black hole mass. But at the same time also many massive 
black holes are excluded that accrete only at low rates. For low 
spheroid masses only positive outliers from the relation can be 
detected, while for increasing spheroid masses also black holes 
on and below the relation can still be detected. This causes the 
bias at low MB„i ge . 

In a realistic AGN sample there may be additional selection 
effects influencing the distribution of MB U i ge . For example, usu- 
ally a minimum host-to-nucleus ratio is required for the AGN 
host galaxy to be detected. This will correspond to some sort 
of cutoff in Mfiuige, excluding low mass hosts. This additional 
selection effect will decrease the average bias for the sample, 
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Fig. 5. M. -Msuige diagram for a Monte Carlo simulation of a 
local AGN sample. The differently sized and coloured symbols 
correspond to different luminosity limits applied to the sample. 
Blue, green, red, and yellow symbols show the sample when 
culled at a bolometric luminosity of / m ; n = 45, 45.5, 46, and 
46.5 in logarithmic units. T he solid black line sh ows the input 
relation for the sample from lHaring & Rixl (120041) . 

as the excluded objects are the most heavily biased. However, 
this will only be of importance if objects are actively eliminated 
during the sample construction based on these considerations. 
Therefore, it is important to understand and quantify the selec- 
tion effects as precisely as possible for a realistic prediction of 
the average bias of the given sample. 

3.2. Flux limited samples 

So far, we have only discussed purely luminosity limited sam- 
ples. This approximation is valid for samples spanning small 
ranges in redshift. However, more generally we have to consider 
samples defined by limitations in apparent flux. (Note that lumi- 
nosity limits may be present in addition, because the AGN need 
to be identified as such and separated from their host galaxies.) 
We first define the multivariate distribution function ^(s, p., I, z), 
as 

¥ ( J, n, /, Z) = Oil, Z) PacO", z) pxQ-H, Z) gQl I s, z) ®s(s, z) ^ .(21) 

dz 

The bivariate distribution function for a flux limited sample is 
then 

¥ (s,/i) = ^ r jj V (s,(i,l,z)dzdl (22) 

where AV C is the comoving volume within the redshift range. 

At this place we ignore possible redshift dependences of the 
individual distribution functions, postponing their discussion to 
section l4~2l We also assume that the selection function is purely 
defined by an AGN luminosity limit and an AGN flux limit, i.e. 
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we neglect any potential dependence of the selection function on 
the bulge property s. The selection function Q is then given by 







for / > & / > / m 



else 



(23) 



I 



log(47rt/ / 2 ) in logarithmic 



where / is the bolometric flux / 
units. 

At high redshifts and for a narrow z range, this is almost 
identical to the luminosity limited case. However, it describes 
the more general case for realistic observations and remains 
valid when the sample covers a wider z range, including the 
low z regime. The application of this form to low redshift AGN 
samples is illustrated in the next subsection. Consequences for 
higher z are discussed in section|4] 

3.3. Application to the reverberation mapping AGN sample 

We have demonstrated above that already the low z AGN popu- 
lation will be affected by the details of the sample selection. The 
main difference compared to high z is that for AGN in the local 
universe the luminosity limit is generally low and therefore the 
bias is less severe. Nevertheless, it is present, and as a test case 
we briefly estimate its influence on the AGN sample with black 
hole masses determined via reverberation mapping. 

The technique of reverberat i on m apping (e.g. 
iBlandford & McKed 1 19821 iPetersonl 1 1993b provides the 
most precise black hole mass estimates for type 1 AGN. It 
builds the foundation for the virial method through establishing 
a scaling relation between the size of t he broad line region and 
the con tinuum luminosity dKaspi et al.Ll200d.l2005tlBentz et all 
l2009al) . However, this method currently only provides mea- 
surements of M. up to a scale factor depending on the (largely 
unknown) geometry and dynamics of the broad line region. 
The usual approach to determine this scale factor and thus 
fix the virial mass scale is to scale the reverberation mapped 
black hole masses of galaxies with stellar velocity dispersion 
measurements to the local M.- < t, rela tion of quiescent galaxies 
(lOnken et all l2004t IWoo et all 120101) . However, this implicitly 
assumes that the reverberation mapping (RM) sample follows 
the same relation as quiescent galaxies. If its M.-cr t relation is 
biased by selection effects then this bias will propagate into the 
absolute calibration of the virial mass estimates. 

Properly defining the selection function for the RM sample 
is difficult, due to the heterogeneous selection of the objects. 
The RM sample as such is certainly not a well defined sam- 
ple. Nevertheless, in order to estimate potential systematics that 
may be inherent in the sample, we tried to estimate the selection 
function, as a kind of first order approximation. The full sam- 
ple covers a wide range in luminosity, from bright quasars to 
moderate-luminosity Seyfert galaxies. Especially for the Seyfert 
galaxies there is no clear AGN luminosity dependence for the 
selection. However, there must be an implicit limit on the AGN 
luminosity and also on the AGN flux. First of all, the AGN has to 
be luminous enough to classify the galaxy as harbouring an ac- 
tive nucleus. Furthermore, the continuum and broad lines have 
to be sufficiently bright to enable spectroscopy and the measure- 
ment of a reliable time lag. For simplicity, we relate the low lu- 
minosity limit to the faintest AGN in the sample of Bentz et al. 
J2009ah . NGC 4051, with log ^agn^ 100 A) = 41.9. This ap- 
proximately corresponds to logLboi = 42.8. Equally, we use the 
object with the lowest flux, PG 1229+204, to define the flux 
limit. This yields a bolometric flux limit of / m i n = -10.4 in log- 
arithmic units, which is also in good agreement with the mag- 




Fig. 6. Predicted M, - MB U i ge diagram and its projections for a 
local type 1 AGN sample, simulating the reverberation mapping 
sample. A flux limit in the redshift interval (0.003, 0.15) and a lu- 
minosity limit are applied to the sample, based on the observed 
range for the reverberation mapping AGN. For this sample a 
mild bias of 0.09 dex is estimated. 



nitude limit of the Palomar-Green (PG) qua sars in the sample 
(ISchmidt & Greenl 119831: iKaspi et all l2000h . We adopt a red- 
shift range of 0.003 < z < 0.15, which contains the vast majority 
of the RM AGN. We also tested the effect of extending the range 
to z < 0.3, including the full RM sample. 

If we assume only a luminosity limit, the total sample bias 
(Ajj.) predicted by our model would be negligible (see Fig 0). 
Indeed, when omitting the flux limit we find (Afi) = 0.003, as a 
lower limit to the bias. However, the luminosity limit is only im- 
portant for z < 0.008. At higher z the flux limit dominates the se- 
lection function. Incorporating the flux limit we get (A/j) = 0.09 
for z < 0.15 ((Ajj.) = 0.12 for z < 0.3). If we relax our conserva- 
tive flux limit, the expected sample bias would decrease slightly. 
Decreasing the assumed luminosity limit has almost no effect, 
as we are dominated by the flux limit. Increasing the luminosity 
limit will slightly increase the bias. 

In Fig. [6] we show the predicted M,-M^ u \ gf , diagram. At the 
high mass end no bias is expected, whereas a mild bias at low 
masses is present. Thus, the M,-M^ u \ ge relation is slightly af- 
fected. 

We further tested this result with Monte Carlo simulations. 
We restricted our Monte Carlo sample to the same luminosity 
limit and flux limit and fitte d the culled sample w ith a maximum 
likelihood method (see e.g. iGiiltekin et all 120091) . with slope b, 
normalisation a and intrinsic scatter <x as free parameters. From 
the restricted sample, we fitted 1000 random subsamples of 100 
objects each. The distribution of the free parameters are shown 
in Fig. [7] We found mean values of a — 8.31, b = 1.12 and 
<x = 0.29, compared to input values of a — 8.2, b = 1.12 and 
cr = 0.3, i.e. we recover the predicted mean offset and confirm 
that the slope and scatter are not strongly affected. 
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Fig. 7. Results of Monte Carlo simulations of the reverberation mapping AGN sample. Panels from left to right: (1) M.- 
Mfluige relation for a Monte Carlo realisation of 100 objects. The dashed line shows the input relation, the solid line the mean 
best fit for 1000 realisations. (2) Distribution of recovered zero points of the M.-Ms u i ge relation. The horizontal line shows the input. 
(3) same for the slope of the relation. (4) same for the intrinsic scatter in the relation. There is an offset in the zero point of ~ 0.1 dex. 
The slope and intrinsic scatter are not strongly affected. 



Therefore, based on our model assumptions, we estimate 
that the reverberation mapping sample may be biased towards 
a high M./Mflutge ra tio by ~ 0.1 dex. This would correspond to 
an underestimation of the virial method by the same amount, 
when usin g the scale factor normalised to the quiescent M.- 
o% relation (lOnken et al.Ll2004tlwbo et all 1201 Oh . However, due 
to the poorly defined selection function of the RM sample this 
can only be taken as a qualitative evaluation of possible selection 
effects. 



3.4. On a black hole mass uncertainty bias 



IShen&Kellvl(l2010l) discussed an additional bias that will affect 
AGN samples with virial mass estimates. This is a Malmquist 
type bias that is expected to arise from the steep turn-down of 
the active BHMF and the intrinsic uncertainty of virial BH mass 
estimates. They argu e that this i s indep endent of the luminos- 
ity bias discussed bv lLauer et al.l d2007l) and thus adds an addi- 
tional bias to the observations. In our framework we are directly 
sampling from the active BHMF, thus we are already taking this 
effect into account. If we introduce an uncertainty in the mass 
estimates this will only increase the scatter, but it will not af- 
fect the mean relation or (Afi). This can again be illus- 
trated by the Monte Carlo simulations. Adding a scatter in the 
BH masses does not affect our sample selection, as it is by con- 
struction based on luminosities, which can be measured with rel- 
atively high precision. Thus the uncertainty just moves the black 
hole masses symmetrically towards lower or higher values, with- 
out affecting the mean. 

It is also straightforward to see this from the distribution 
function y V (s,fi). We assume that the virial mass estimate (the 
observed mass /j. q ) is given by a log-normal probability distribu- 
tion with mean fi and dispersion <x vm , the uncertainty in the virial 
mass estimate, 



g(p IaO = -£= exp-j — 



V2tto 



2oi 



(24) 



The bivariate distribution function for bulge property and virial 
black hole mass is then 



The mean relation is 

, ws rr^og(/"ol/u)^o(*,^)dyud// J fi y ¥ (s,ij)diJ. 

</U(s) = :H 77 = — r , (26) 

JJ g(Mo I /") ^oCs.W dMu J dyu 

when integrating over /j Q . This is identical to Equation ©, so no 
additional bias is introduced by the virial mass uncertainty. 

The conceptual diffe rence between our results and the work 
bv lShen & Kelivl (1201 Ol) is, firstly, that we do not treat the black 
hole mass bias independently, but within the total bias budget. 
Secondly, we account for the fact that in observations we do not 
sample directly from the true BHMF, but only from the apparent 
BHMF affected by the relevant luminosity limits. As discussed 
by SW10 in detail, this directly observable distribution suffers 
from incompleteness at the low mass end and therefore turns 
over towards low M.. 

In Fig. [8] we show the result for a luminosity limit of / m i n = 
47, assuming <x vm = 0.3 dex. While the black hole mass distri- 
bution is broadened, the mean relation is unchanged. 



3.5. On a spheroid uncertainty bias 

At least at high redshifts, not only the estimated black 
hole masses are significantly uncertain, but also the observed 
spheroid properties. For example, the emission line width of 
[Om] has been used several times as surrogate of cr». However, 
this is at best on a verage a reliabl e estimator, with a dispersion of 
at least ~ 0.2 dex jNelsonl [20001) . Furthermore, spheroid masses 
estimated for AGN host galaxies suffer from several uncertain- 
ties, including the difficulties of performing a proper AGN-host 
decomposition, of deblending bulge and disk components, and 
of converting host galaxy luminosities into stellar masses, using 
either a fixed mass-to-light ratio or colour information. These 
uncertainties may well accumulate to a total error of similar 
magnitude as for the black hole masses, or more. In that case, 
again a Malmquist type bias is possible. 

We model the observed spheroid property s Q by a log-normal 
probability distribution with mean s and dispersion o~ s , analog to 
Equation d24l i. The bivariate distribution function is then 



T (j,Mo)= sGUol/W^d/i. 



(25) *F (w) 



g(s \s) x P (s,fi)ds, 



(27) 
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Fig. 8. Predicted M. - MB U i ge diagram and its projections for a 
type 1 AGN sample, with luminosity limit of l m \ n = 47 (as the 
right panel of Fig. [3]). For the black hole masses an uncertainty 
of 0.3 dex is assumed, representing the intrinsic uncertainty in 
the virial method. A thick solid line shows the observed M.- 
A^Buige relation for the observed M., which is identical to the re- 
lation without uncertainty in the virial mass. Only the scatter in 
the relation is increased. The red dashed dotted line shows the 
mean relation (s)(/i ), while the dotted red line gives the mean 
relation (s)(ji), i.e. without measurement uncertainty in the black 
hole mass. 



and the mean relation is 

fffig(s \s) *P (i, fi) ds d/i 



ffg(s \s)W (s,fi)dsdfi 



(28) 



The total sample bias (A/j) is not affected by the uncertainty, as 
can easily be verified. However, the predicted M.-bulge relation 
is affected, as shown in Fig. [9] At the low mass end the mean 
M. increases compared to the prediction without scatter, while 
at the high mass end it decreases. To understand this trend we 
have to recall that we are sampling from the luminosity-limited 
AGN host galaxy spheroid mass function. Due to the luminosity 
limit the density of sampled objects decreases at the low mass 
end, similarly as in the active BHMF (cf. SW10). Because of 
this decrease at the low mass end more objects are scattered 
from slightly larger masses to lower masses than the other way 
around. This produces an excess of higher true s at a given ob- 
served s , which in turn leads to an increase in the average y. at 
the given s a . At the high mass end the reverse happens. Due to 
the steep decrease towards higher masses, more black holes are 
scattered to higher s a , having on average a lower ju, as they in- 
trinsically have a lower s than observed. In total, the observed 
M.-bulge relation flattens, without changing its normalisation. 

For the mean spheroid mass at a given black hole mass the 
opposite effect happens. Uncertainties in the spheroid masses do 
not change the mean relation (s )(/i), but uncertainties in the 
black hole masses can strongly affect the mean relation {s)(/u ), 
causing a steepening of the observed relation. This is indicated 
by the red lines in Figures [8] and [9] 
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Fig. 9. Same as Fig. [8] but with measurement uncertainty of 
0.2 dex in the spheroid mass and no uncertainty in the black hole 
mass. The thick solid line shows the M.-MB U i ge relation for the 
observed MB U i ge values, while the dotted line shows the same re- 
lation but for the true MB U i ge , i-e. for the case of no measurement 
uncertainty. The sample bias is not affected, but the slope of the 
M.-MBui e e relation flattens due to uncertainties in MB U i ge - The 
red dashed-dotted line shows again the mean relation (s )(ju), 
which is not affected by uncertainties in MB U i ge - 



Finally, we note that in general the convolution with the un- 
certainty should really be done with ^(i,^) rather than with 
TVs,^), i.e. before applying the selection criteria. In our case of 
a purely luminosity limited sample, the selection function does 
not depend on s a or fi . Therefore, our approach is justified and 
no additional bias is introduced for the sample in total. In case of 
a more complicated selection function that depends on s or yu 
this might no longer be the case. An additional bias by the mea- 
surement uncertainty then becomes possible which will depend 
on the details of the selection function. 



3.6. The slope in low-redshift AGN samples 

The black hole - bulge relations determined for low-redshift 
AGN samples usually plays an important role in the search for 
possible evolution in the relations. But they are also interesting 
in themselves, as it is not at all clear from a physical perspective 
that AGN host galaxies should follow exactly the same relation 
as largely inactive galaxies. Notice that the current calibration 
of the virial method ensures by design that AGN and non-AGN 
relations have no global offset from each other. But it has been 
noted that the slope of the M.-bulge relation for some local AGN 
samples appeared to be flatter than in the relation of black holes 
with dynamical black hole mass measurements dGreene &HoL 
[20061: iKmTetail l200l iBentz et all 120095 IWoo et all 120101) 
We were interested to see whether this flattening can be ex- 
plained within our framework by selection effects, or whether an 
intrinsic difference in the relations between active and inactive 
galaxies has to be invoked. 
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Fig. 10. Comparison of our bivariate probability distribution ^(s,^) (contours), assuming reasonable numbers for the selection 
function, with observational data of low-redshift AGN. The thick black line shows the mean relation (iu)(s), the red dashed dot- 
ted line shows {s)(/j), the black dashed line gives the true relation. The dotted lines shows the mean relation (ij)(s) including 
measurement uncertainties in s of 0.05 dex for cr* and 0.1 dex for Ly. L eft panel: M.-cr, r elation. We show stellar velocity disper- 
sion measurements for the reverberation mapping sample (blue circles, IWoo et all 1201 Oh . and a local SPSS sample (green stars, 
iGreene & Holl200"6h . Right panel: M. - Ly rel ation for the rever beration mapping sample (blue circles, iBentz et all l2009bh and a 
local QSO host galaxy sample (green triangles. iKim et al 1 120081) . 



In Fig [10] we show the AGN samples for which such evi- 
dence for a flatter slope of the M. -bulge relations has been di- 
agnosed. We compare them with predictions of the bivariate dis- 
tribution function *T (s,//), adopting plausible numbers for the 
selection function. For the M.-c r„ relation we used t he stellar ve- 
locity dispersion function from ISheth et all (120031) and assume 
our approximation of the selection function for the reverbera- 
tion mapping sample from section 1331 For the M, — Ly relation 
we ado pted the galaxy luminosity function from lKochanek et all 
(1200 11) . converted to t he V-band and into a spheroid luminos- 
ity function following iMarconi et al. 1 J2001 . To roughly model 
the selection function we assume an AGN luminosity limit of 
logLboi = 10 45 ergs s -1 . 

Adopting these selection functions, we show in Fig Q7J the 
predicted mean relations (fi)(s) as thick solid lines. The bias 
caused by the AGN selection depends strongly on s and is largest 
at low s. Thus the relation does not only change in offset but 
also flattens compared to the true relation. Furthermore, mea- 
surement uncertainties in the bulge property lead to a further 
flattening of the slope, as discussed in 13.51 Indeed, for active 
galaxies the bulge property is usually measured with larger un- 
certainty than for inactive galaxies. We computed the expected 
mean relations (ij)(s) adding reasonable values for the measure- 
ment errors in the respective spheroid properties to the bivariate 
distributions. These are shown in FigQJJas dotted black lines. 

Considering our rough approximations to the respective sam- 
ple selection functions we find a very good agreement between 
our prediction of a flatter M.-bulge relation and the observations. 
Thus we argue that these selection effects are able to reproduce 
the observed flattening in the slope. There is no reason to assume 
that AGN have generally a different black hole - bulge relation 
than inactive galaxies. 

We note in passing that the M.-Ly relation seems to disagree 
with our bivariate probability distribution at the high mass end. 
This may reveal a problem with our assumption, but it may also 
just indicate that the M.-Ly relation is not the same for inac- 
tive galaxies and AGN host galaxies. In fact that is even ex- 
pected to be so, since it is known that host galaxies of high lu- 
minosity AGN have on average younger stellar populations (e.g. 



Kauffmann etall l2003t Uahnke et all BOOl IVanden Berk et all 
20061) ' compared to inactive galaxies. Thus, if they obey the same 
M.-Mbuige relation, they actually should deviate in the M.-Lfiuiee 
relation. 



4. Evolution in the M.-M Bu i g e relation 
4. 1 . Evolution in a flux limited sample 

How will selection effects influence studies that test for redshift 
evolution in the black hole-bulge relations? If exactly the same 
sample selection criteria are applied at high z as for a z ~ com- 
parison sample, no bias will be present, at least to first order (see 
section l4~2l i. Any change in the M.-bulge relations could then be 
taken as evidence for their evolution. However, this simple situ- 
ation is usually not given. The local comparison is commonly 
provided either directly by the quiescent relation, or by a lo- 
cal type 1 AGN comparison sample. While the former clearly 
possesses different selection criteria, also the latter is not auto- 
matically selected in the same manner. This has to be taken into 
account when comparing low and high redshift results. 

To illustrate this point we derive the sample bias for a flux 
limited sample at a given redshift, 



<4p>(z) 



ffj(pt-a-bs) <F (s, fx, I, z) dl dfi ds 
Jff y ¥ (s,fx,l,z)dld/xds 



(29) 



with the multivariate distribution function ^(s,/^, l,z) given by 
Equation (|2TT >. We thus assume these distributions, for the mo- 
ment, to be non-evolving with redshift (we discuss the more gen- 
eral case below). Equation ( f29b is equivalent to the bias for a 
luminosity limited sample, but with a redshift dependent lumi- 
nosity limit. 

Instead of using a bolometric flux limit as before, we now 
adopt an optical flux limit, for example in the B band. For sim- 
plicity we assume a simple power law K-corr ection with spectral 
index a = -0.44 (IVanden Berk et all l200ll) . The B band lumi- 
nosity is converte d to bolometric lumin osity using the bolomet- 
ric corrections of IMarconi et all (120041) . In Fig. [TTJwe provide 
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Fig. 11. Redshift evolution of the selection bias for a magnitude 
limited AGN sample that follows the local distribution functions 
at all redshifts. The dashed, dotted, solid, and dashed dotted lines 
show B band magnitude limits of 19, 20, 22, and 24 mag, respec- 
tively. 

the expected sample bias at a given z for different apparent mag- 
nitude limits. The deeper the respective survey, the lower is the 
expected bias. For a fixed flux limit, the expected bias increases 
with redshift, which can mimic an evolutionary trend. 

Furthermore, to have a comparable bias between a local sam- 
ple with, say, B^ ss 19 mag and a sample at z » 2, the high 
redshift sample has to be complete down to B m j n « 24 mag. 
The use of a local AGN sample as comparison sample for higher 
z studies is already an improvement over the use of the quies- 
cent black hole sample, but it does not ensure the absence of a 
selection bias. Ideally, the local comparison sample should be 
matched in AGN luminosity to the high z sample. 

On the other hand, if the selection function of the respective 
sample at high z is known, the expected bias can be computed, 
and any clear offset from this prediction could be interpreted as 
evidence for evolution in the M. -bulge relations. However, we 
show in the following subsection that unfortunately even this ap- 
proach is not fully unbiased. 

4.2. AGN evolution biases 

So far we ignored any possible explicit redshift dependence of 
the underlying distribution functions in Equation (OTT i. However, 
at least some of them must evolve with redshift. The stellar mass 
function, and thus also the spheroid mass function, are certain!) 



changing with z (e . g^ Bundvetafl 120051: Franceschini et 



inly 

"all 



l200ia iPozzettietall l2007t lllbert et all 120101) . Furthermore, it 
is well known that the AGN population itself strongly evolves 
apparent in the evolution of the AGN luminosity function (e.g 



Ueda et all l200l[ lHasinger et all l2005t iRichards etall 12006b 
Bongiorno et all 20071: ICroom et all 120091) . Not only the nor- 
malisation and typical luminosity change with z but also the 
shape of the luminosity function. At low z the faint end of 
the QSO luminosity function steepens, known as 'AGN cos- 
mic downsizing". While in the local uni verse the QSO lumin os- 
ity function shows only mild curvature (I Schulze et all l2009h . a 
prominent break is present in the high z luminosity function. 
This directly implies evolution in the active black hole mass 



function, the Eddington ratio distribution function, or most prob- 
ably in both, leading to a change of the predicted bias with red- 
shift, even for a fixed luminosity limit. Practically speaking it is 
generally not straightforward to judge if an observed evolution 
in the M, -bulge relations is caused by evolution in the intrinsic 
relations, or whether the trend is spurious and caused by evolu- 
tion in the underlying distribution functions. 

A prediction of the sample bias at given redshift requires, 
apart from a well defined selection function, also knowledge 
about these underlying distribution functions, the spheroid mass 
function, active BHMF and ERDF In the local universe these 
distributions are at least reasonably well established (SW10), 
but for the high redshift universe the situation is currently 
less comfortable. The distribution function of velocity dis- 
persions is reasona bly well established for the local u niverse 
dSheth et all l2003[). but rat her uncertain beyond that dChaeL 
120101 : iBezanson et all 1201 ll). The total stella r mass function is 
observed up to z ~ 4 ( Fontana et all 120061) . but galaxy mass 
functions for different morphological types have been deter- 
mined only up to z ~ 1.4 dBundv et alll2005l : lFranceschini et all 
120061: lllbert et all 1201 Oh . The spheroid mass function itself is 
again basically unknown for higher z. To enable an analysis of 
how evolution influences the biases, we now derive very rough 
estimates of these distribution functions for higher redshifts. 
These evolving distributions are meant as illustrative only; a 
more detailed investigation would be beyond the scope of the 
current paper. 

An upper limit to the spheroid mass function is given by the 
total stellar mass fun ction. We use the para metric fit to the stel- 
lar mass function by iFontana et all d2006h for this purpose. A 
lower limit is given by the mass function of elliptical galax- 
ies. To derive this distribution for arbitrary z we assume the 
same ellip t ical-to -total ratio as for the local mass function from 
iBell et al] d2003h . and apply this correction to the total mass 
function at higher redshift. This is clearly an oversimplification, 
as the relative contribution of elliptical galaxies seems to de- 
creas e with increasing z, at least to z ~ 1 (e.g. iBundv et all 
120051) . However, at higher redshifts the low mass end that is 
most affected by this correction is not well determined and may 
be underestimated if also the distribution of mass-to-light ratios 
changes. For the present purpose we hold that our simple approx- 
imation serves as a reasonable lower limit of the spheroid mass 
function. Additionally, we define a case in between these lim- 
its, using the same conversion of total mass function to spheroid 
mass function as assumed above for redshift zero (section [2~TI ). 

There has been significant progress in the determination 
of the active BHMF and Eddington ratio distribution func- 
tion at high redshifts i n the last years (Vestergaard et al., 2008; 
Vestergaard & Osmen. l2009t iKellv et all 120101: IShen & Kelly 



201 lh . But these results mainly cover the bright end of the lumi- 



nosity function and thus the high mass end of the BHMF, while 
the low mass end is still poorly determined; also the systematics 
are not fully understood. Here we restrict ourselves to the use of 
a mass function that is consistent with current observations. 

We now used the local active BHMF and ERDF as well as 
the redshift evolution of the type 1 AGN luminosity function 
(LF) as constraints to predict the distribution functions at higher 
redshift. To achieve a smooth redshift evolution of the BHMF 
and ERDF we explored several arbitrary, but reasonable evolu- 
tion models for these distribution functions. We fixed the zero 
point of the two distribution functions to their local values (as 
given by SW10), and fitted their redshift evolution to the opti- 
cal LF (employing Equation [T2b. For the LF we used the lumi- 
nosity dependent density evolution model from lBongiorno et all 
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Fig. 12. Reconstructed AGN distribution functions at three representative redshift derived from fitting the observed AGN luminosity 
function to a redshift evolution model for the BHMF and ERDF. For the BHMF a mass dependent density evolution model is 
used, and no evolution in the ERDF is assumed. Left panel: Best fit active black hole mass function. Middle panel: Eddington ratio 
distribution function. The shape has been fixed, the normalisation is determined by the space density of the black hole mass func tion. 
Left panel: AGN luminosity function. The red lines show the type 1 AGN luminosity functions from lBongiorno et al.l d2007l) . the 
black lines are our best fit to them. 



(I2007h . which is based on the SDSS at the bright end and on the 
VVDS at the faint end. Optical B band magnitudes were con- 
verte d to bolometric lumin osity using the bolometric corrections 
from lMarconi et al"1(l2004h . Alternatively we also explored mod- 
elling the BHMF and ERDF directly from the observed LF in in- 
dividual redshift bins, which lead to the same qualitative results. 

As long as only the LF is used as constraint there is clearly 
a degeneracy between evolution in the BHMF and in the ERDF, 
respectively, while presumably both are evolving. We tried to 
bracket this degeneracy by exploring two extreme cases. Firstly, 
we assumed a constant ERDF and let the BHMF change with z. 
Secondly, we fixed the active BHMF and let the ERDF evolve 
with z- 

The first case, a non-evolving ERDF, can be regarded as an 
upper limit to the expected AGN evolution bias. The downsizing 
of the AGN LF directly corresponds to a downsizing in the ac- 
tive BHMF. In Fig.[l2]we show the best fit BHMFs, ERDFs, and 
the corresponding AGN LFs, for a few representative redshifts. 
Here we adopted a mass-dependent density evolution model to 
model the redshift depende nce of the active BHMF , inspired by 
the LDDE model used by iBongiorno et al.l d2007l) to parame- 
terise the evolution in the AGN LF. In order to reproduce the 
observed AGN LF, the space density at the high-mass end has 
to increase substantially. Thus, the mass dependence of the ac- 
tive fraction is reduced, which in turn changes the active fraction 
bias. Because this effect tends to counteract the AGN luminosity 
bias, a reduced active fraction bias leads to an increase of the 
expected overall sample bias compared to the redshift zero case. 
This is illustrated in Fig. Qj] where we show the expected sam- 
ple bias (A/j) as a function of the applied lower luminosity limit 
for various redshifts. 

The second case, a fixed BHMF and an evolving ERDF, 
provides a lower limit to the expected AGN evolution bias. In 
this scenario the active fraction is only mildly evolving with 
z, through the evolution of the spheroid mass function. On the 
other hand, the ERDF needs to evolve strongly, in the sense of 
an increasing average log A , to satisfy the constraints from the 
AGN LF. For a higher average log A the same luminosity limit 
corresponds to a lower average black hole mass limit, thus re- 
ducing the expected bias. Therefore, this assumption will pre- 
dict the lowest sample bias. However, the best fit model en- 
forces an unreasonably high space density of objects accret- 
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Fig. 13. Expected sample bias as a function of the lower lumi- 
nosity limit for the same representative redshifts as in Fig. Q~2] 
The dotted line shows the z ~ prediction, the other lines show 
predictions at higher z- The shaded area incorporates the uncer- 
tainty on the spheroid mass function, by using the total mass 
function as upper limit and an estimate of the elliptical mass 
function as lower limit. 



ing at super Eddington rates and also provides a poorer fit to 
the AGN LF. Furthermore, it dis agrees with the pictu r e of anti- 
hierarchical black hole growth dMarconi et all 120041 iMerloniL 
120041 : iMerloni & Heinzl 120081) . Thus, this case is not physi- 
cally plausible. We also experimented with intermediate evolu- 
tion scenarios, where a mild evolution of the ERDF balances 
the increasing bias through the evolution of the BHMF to some 
degree. For all physically plausible scenarios, we found at least 
some increase of the bias compared to the local case. We con- 
clude that the estimated sample bias for the local universe is a 
lower limit for high z samples. 
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Fig. 14. Redshift evolution of the selection bias for a magnitude 
limited AGN sample, assuming our estimate of the z evolution 
of the underlying distribution functions, as discussed in the text. 
The shaded regions incorporate the uncertainty on the spheroid 
mass function. 



In Fig. [14] we show the expected bias at a fixed magnitude 
limit as a function of redshift, equivalent to Fig. QT] but now 
allowing for an evolving active BHMF and a constant ERDF (our 
case 1 above). As discussed, the resulting bias can be taken as 
an upper limit. Within the assumptions made in our model, only 
observations that show a clear excess on top of this predicted 
bias would constitute firm evidence for real evolution in the M,- 
bulge relations. 

While evolution is usually expressed in terms of the normal- 
isation of the relation , it is also pos s ible th at there is evolution 
in its intrinsic scatter dMerloni et all l20loh . The strength of the 
sample bias is strongly affected by the amount of intrinsic scatter 
in the M.-bulge relations. Therefore, an offset of the M./Mbui ee 
ratio on top of the expected selection bias can equally be inter- 
preted as being caused by an increased bias due to an increased 
intrinsic scatter, rather than as a true offset in the zero point of 
the respective M.-bulge relation. To illustrate this effect, we ar- 
bitrarily assume a redshift evolution for the intrinsic scatter as 
<t(z) = cr z=0 (l + z) 05 with 0"j=o = 0.3, giving approximately 
cr(z = 2) » 0.5. In Fig.[15]we show the resulting redshift evolu- 
tion of the bias for fixed magnitude limits. The predicted scatter 
is strongly enhanced in this case, especially at high redshifts. 

An increase in the intrinsic scatter is also some kind of evo- 
lution in the relations. However, its interpretation is of course 
different from a change in the zero point. Whereas the latter 
would suggest a phase where black holes grow stronger than 
galaxies, or the other way around, the former would be consis- 
tent with a true coeval growth on average, in which the correla- 
tion is tightened from a rather loose one to the tight correlation 
we observe today. Indeed, such an evolution in scat ter towards a 
tighter relation is expected in various model s (e.g. iPend, 120071: 
IVolonteri & Nataraianl 12009b lLamastra etall 1201 Oh . Given suf- 
ficiently high quality data and large enough samples, it should 
be possible in the future to estimate the amount of evolution of 
the intrinsic scatter in the relation. 
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Fig. 15. Same as Fig. [14] but assuming an evolution of the intrin- 
sic scatter in the M.-M Bu i g e relation of <t(z) = 0.3(1 + z)° 5 . 



5. Discussion 

Several observational studies of type 1 AGN samples found ten- 
tative evidence for an increase in the M./Mn u \ ge ratio with red- 
shift. This increase is often parametrised as an evolutionary be- 
haviour M./MRni^ oc (1 +z) y , with values for y up to 2.1 reported 
dMcLure et al.Ll2006l) . As discussed in this paper, these observa- 
tional studies are invariably affected by selection effects. We now 
apply our formalism to a few published studies. We estimate ap- 
proximate distribution functions and predict the resulting sample 
biases. Note that we do not aim at deriving detailed corrections, 
or at presenting an exhaustive discussion of the relevant liter- 
ature. Rather we wish to highlight the magnitudes of possible 
systematic biases. 

Specifically, we explored two scenarios for the distribution 
functions that serve as upper and lower limits for the magnitude 
of our predicted selection bias. 

• Model 1: We used the local distribution functions presented 
in section 13.11 throughout the entire redshift range, i.e. we 
ignore any effect of redshift evolution. As discussed above, 
this will serve as a lower limit to the expected bias even if 
the scenario is unrealistic. 

• Model 2: Here we incorporated redshift evolution in the un- 
derlying distribution functions, in particular in the AGN pop- 
ulation. We used the model discussed in section 1431 assum- 
ing a mass dependent density evolution for the BHMF and 
a non-evolving ERDF. This model provides an approximate 
upper limit to the selection bias. 

5.1. \Merlonietai\l20Tdi) 

A relatively large an d well-defined sample was investigated by 
iMerloni et al.l d2010h . Th ey used an 7 band limited set of objects 
drawn from zCOSMOS (iLillv et all I2007I) in the redshift range 
z = [1.06,2.19], with 7ab < 22.5, for which they estimated stel- 
lar masses by fitting single stellar population models to the ob- 
served multiband SEDs. The results are shown in Fig. [16] While 
the individual data points show considerable scatter, they are 
on average of f set fro m the z = relation by (Ayu) b s = 0.34. 
IMerloni et al.l (l2010h considered also selection biases, largely 
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Fig. 16. Le ft and middle panels : Predicted M,-M% a \ ge probability distributions for the samples of iMerloni et al.l (1201 Ot) (blue filled 
circles) and lJahnke et al.l d2009t) (yellow squares), using our model 1 (left) and model 2 (middle) assumptions about the distribution 
functions. In each panel, the black solid line shows (fj)(s), the red dashed-dotted line shows (s)(ju), and the thin black dashed line 
delineates the reference z ~ relation. Right panel: redshift dependence of the observed M. /MB U i ge ratios, expressed as offset A/i 
from the z — relation, together with the predicted evolution of the sample bias. The symbols are as before, plus the large open 
circles that show the mean offset for the Merloni data in three redshift bins, and the large open square showing the mean of the 
Jahnke data. The dashed line shows our prediction for model 1, the solid line with the shaded area shows our prediction for model 2, 
including the uncertainties on the spheroid mass function. The dashed-dotted line is for model 1, but assuming evolution in the 
M.-M Balge relation M./M Bu i g e « (1 +z). 



following the argument given by (lLauer et al.ll2007|) . as one pos- 
sible origin of this offset. They concluded that their observations 
indicated either evolution in M. /M% u \ ge or a significant increase 
of intrinsic scatter on the M.-M Bu \ ge relation. 

We modelled their sample selection and applied the two 
above mentioned models for the distribution functions to predict 
the expected sample bias. The correspondingly modified bivari- 
ate probability distributions are overplotted in the left and mid- 
dle panels of Fig. [16] Neither reproduces the distribution of data 
points very well, but it seems that an intermediate model might 
actually do a reasonable job. In the right-hand panel of the figure 
we compare the expected average sample offset due to selection 
biases, as a function of redshift, with the data points (individual 
and also binned into three redshift ranges). The binned mean off- 
sets agree quite well with the range of values covered by our bias 
model. The prediction for a M. /MB U i ge ratio evolving with red- 
shift as y — 1 , on the other hand, lies slightly above the binned 
data. 

We pointed out in Sect. 13.11 that the relation of mean stel- 
lar masses at given black hole mass is to first order unaffected 
by the AGN selection bias and depends mainly on the spheroid 
mass function. We binned the data points in M. accordingly 
and compare this binned relation in Fig. [17] with the prediction 
for the null hypothesis (i.e., consistency with the local relation). 
We also computed the same predicted relations for an evolving 
^•/A^Buige ratio (again with y — 1), and also for a higher in- 
trinsic scatter in M, /Mfluige of 0.5 dex (but no evolution). Both 
curves are offset from the null hypothesis towards higher BH 
masses, albeit at slightly different tilts. We performed a simple 
X 2 test to compare the binned data with the predicted distribu- 
tions. The lowest value of x 1 , and in fact the only one close to 
a reduced^ 2 of unity, was obtained for the null hypothesis. We 
thus conclude that this sample provides no statistically signifi- 
cant evidence for evolution i n M. /Mb„i P r or its sc atter, and that 
the apparent offsets found bv lMerloni et al.l ( 120101) are consistent 
with the expected magnitude of AGN sample selection effects. 
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Fig. 17. Comparison of predictions for the relation of mean 
stella r mass at giv e n bla ck hole mass (s)(fj) with the data 
from IMerloni et all d2010l) . The blue circles show their data, 
the red triangles give their mean stellar masses in black hole 
mass bins. These are compared with predictions of (s)(//) using 
our model 2. We employ the case of no evolution in the M.- 
MsuUe relation and an intrinsic scatter of cr — 0.3 (solid line), 
evolution in this relation with gamma = 1 (dashed line) and for 
a higher intrinsic scatter of 0.5 dex but no evolution in the rela- 
tion (dotted line). The black dashed line again indicates the true 
M.-Mfiuige relation. 



5.2. \Jahnke et all 1200$) 

iJahnke et al.l d2009l) studied t he M, - M, re l ation for a sample 
of 10 AGN from COSMOS dHasinger et al.l 120071) . using opti- 
cal and NIR HST images to estimate stellar masses of the host 
galaxies. Their sample is X-ray selected, but the spectroscopic 
follow-up requi rement leads t o a sim ilar optical flux l imit a s 
for the work by IMerloni et afl (1201 Oh (see iTrump et all l2009h . 
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Fig. 18. Observed and predicted redshift evolution of the offset in 
black hole masses for a sample of SPSS quasa rs. The green cir- 
cles show the data from lSalviander et al.l d2007h . using the [O in] 
line width as surrogate of cr», the blue squares show their re- 
sults using [On], and the yellow triangles show the results from 
IShen et all (120081) . The dashed and solid lines shows our model 1 



and model 2, respectively, assuming no intrinsic evolution in the 
M.-cr* relation. The dashed-dotted and dotted lines give the pre- 
dicted sample offsets for an evolution with y — 1 for the two 
models. The shaded areas again incorporate the uncertainty on 
the spheroid mass function. 
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Fig. 19. M. a nd M, measurements of high-luminosity quasars 
at z ~ 3 from I Schramm et all d2008h . The black contours show 
the predicted M,-M Bu \ %e probability distribution using model 2 
for the null hypothesis of no evolution in M. /MB U i ge ■ For com- 
parison, the offset green contours show the expected distribution 
if M. /MBuige were to evolve as strongly as with y = 2. 



As the covered reds hift range is very sim ilar as we ll, we simply 
added the results bv lJahnke etail (120091) to Fig. [16] Uahnke et all 
(120091) found their sample offset Afj. to be close to zero and con- 
cluded that their data are fully consistent with the local relation. 
Fig.[l6]suggests that when accounting for selection effects, these 
data actually might even imply a negative evolution. It is in- 
teresting to note t hat also the small X-ray selected sample by 
ISarria et all (1201 Oh . covering a similar redshift range finds no, or 
even a negative evolution. H owever, notice that both studies (and 
also iMerloni et al.L 1201 0|) estimated only t otal galaxy masses 
instead of bulge masses. Jahnke et alJ {2009) discussed possible 
consequences of this important difference for the true evolution 
of the M.-M Bu i ge relation. 

5.3. \Salviander efa/1 l20ul) 



We now consider the study by ISalviander et al.1 (120071) . who 
studied the M.-cr t relation in the redshift range < z < 1.2, 
using the [Om] line width as surrogate of cr, for z < 0.8, and 
the [On] line width for 0.4 < z < 1.2. Their sample is drawn 
from the Sloan Digital Sky Survey Data Release 3 (SDSS DR3; 
lAbazaiian et alll2005|). which does not constitute a well-defined 
sample (see e.g. iRichards et all 120061) . Furthermore, they ex- 
cluded a large fraction of objects based on quality cuts. To model 
their sample we simply ad opted the flux limit of the SDSS main 
quasar sample, i < 19.1 dRichards et al.L 120021) . In Fig [T8l we 
show our predicted redshift evolution for the sample bias and 
compare them to the obser vations (which we augmented by sim- 
ilar low-z measurements by [Shen_eTaU (|2008|)). For z < 0.6, their 
data points are slightly below, for z > 0.6 slightly above our pre- 
diction for Ayu(z) under the null hypothesis of no evolution in the 



M,-cr, relation. We also computed A/i(z) again for the case of 
evolution with y — 1 . 

In their dis c ussion of possible biases in their results, 
ISalviander et al.l d2007l) attributed ~ 0.25 dex of their offset 
to selection effects, which roughly agrees with our estimate at 
z ~ 0.7. They also discussed another bias caused by their signal- 
to-noise requirement for the narrow lines, which may accumu- 
late to an additional ~ 0.15 dex at high z. Taking this additional 
effect into account, the data appear again to be largely consistent 
with the null hypothesis of no evolution. 

5.4. High redshift luminous QSOs 

The largest apparent offsets from the local M, -bulge relations 
were so far observed for very luminous QSOs at high redshifts , 
z > 3 up to z ~ 6 (e.g. IWalter et all 12004 IWang et all 120101: 
iTargett et all 1201 it) . Because of the limitation to very luminous 
objects, significant selection effects are expected. At these high 
redshifts our knowledge of the required distribution functions 
becomes extremely vague, and we can do no more than guess 
the magnitude of the sample biases. It is nevertheless instruc- 
tive to see how large they can become. We illustrate this by 
using the s mall sample of t hree h yperluminous QSOs around 
z ~ 3 from ISchramm et al. J 2008b. These quasa rs originate in 



the Hamburg/ESO Survey dWisotzki etall l2000t) which has an 
approximate magnitude limit of B < 17.5. Stellar masses were 
estimated from NIR imaging that provided rest-frame B - V 

col ours. 

ISchramm et alJ d2008l) obtained M. /M* values that were an 
order of magnitude larger than found in the local universe. In 
Fig. [19] we superimpose this small sample over the predicted bi- 
variate probability distribution assuming no evolution. The data 
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are well captured by this distribution; in fact, it is expected that 
luminous QSOs would nearly always be found well above the 
local relation. Thus, there is again no evidence even at z ~ 3 that 
the M. -A/Buige relation is any different from today. 

However, the error bars of these and similar measurements 
are considerable, and some level of evolution would clearly be 
equally consistent with the data. To see what the expected signa- 
ture for very strong evolution might be we show in Fig. [T9l also 
the probability distribution for evolution with y — 2, a value that 
has been suggested by some authors in the recent past. In that 
case, offset s of almost 2 order of m agnitude would be expected; 
the data by ISchramm et al.l (|2008) are clearly inconsistent with 
such a scenario. 

The general features discussed with this example are appli- 
cable to other studies of bright QS Os at high z. In p articular, the 
z ~ 6 SDSS main quasar sample dFan et all 1200 ll) has a simi- 
lar luminosity limit, and sample biases at least as strong can be 
expected. However, we make no attempt to guess the relevant 
AGN distribution functions at z ~ 6 and thus have to refrain 
from making a speculative statement about possible biases for 
these highest redshift objects. 



6. Conclusions 

We investigated the ramifications of sample selection effects for 
the observed black hole - bulge relations. Our starting point is the 
bivariate probability distribution of galaxy bulge properties and 
black hole masses, which ideally describes the true underlying 
relation. However, this relation is recoverable only from a strictly 
volume-limited sample with completely measured galaxy prop- 
erties and black hole masses. In real observations the probability 
distribution is inevitably modified. We incorporate the modifi- 
cation of the bivariate distribution by the use of an appropriate 
selection function which represents the observational process of 
sample definition and measurement. In the presence of selection 
effects and scatter, the true M. -bulge relations are generally not 
recovered by the observations. 

Others have recognised the importance of selection effects 
before, and someti mes tried to include them in their bud- 
2007t iMerloni et all 1201 Ot iBennert et all 



get flS alviande r et al. 
120101 Lamastra et al 



20101) . However, different papers focussed 



on different effects, and sometimes with quite simplifying as- 
sumptions. In this paper we go one step further and provide a 
common framework in which all kinds of selection effects on 
the M.-bulge relations can be investigated. We also identified 
some important sources of bias that were not discussed in the 
literature so far. 

Already the widely used M.-bulge relations derived from ob- 
servations of nearby galaxies with dynamical black hole mass 
measurements do not strictly reflect the underlying intrinsic re- 
lations. We showed how the bivariate probability distribution 
is changed due to explicit or implicit selection against galaxies 
where the central black hole sphere of influence is not well re- 
solved. The resulting bias can be quite strong, almost of the order 
of the intrinsic scatter, when a high threshold is adopted. It be- 
comes weaker when the threshold is lowered and more galaxies 
are included, but this comes at the penalty of increasing the ob- 
servational errors. Measuring BH masses in galaxies with poorly 
resolved spheres of influence also puts higher demands on the 
dynamical modelling, such as the inclusion of Dark Matter ha- 
los. 

For the main part of this paper we investigated in detail a 
number of selection effects for type 1 AGN samples, as such 



datasets are essential as probes for redshift evolution in the M.- 
bulge relations. We identified a variety of effects: 

• Active fraction bias: The probability for a black hole to be 
in an active (type 1 AGN) state depends on the black hole 
mass. At least in the local universe, very high mass SMBHs 
are much more commonly quiescent than lower mass black 
holes, which have a higher probability to be actively accret- 
ing. In comparison with an ideal volume-limited set of galax- 
ies, an AGN sample will therefore contain a higher fraction 
of low mass black holes. We have shown that this selection 
effect implies a bias towards a low M./M% u i ge ratio. This 
source of bias has not been discussed before. 

• Luminosity bias: The selection of AGN based on their lumi- 
nosities implies a preferential selection of black holes with 
higher than average masses, especially when the luminosity 
limit is a high one. This effect has already been discussed and 
analy sed by others (e.g. ISalviander et"aU 120071 lLauer et all 
120071) . Our analysis revealed that the luminosity bias de- 
pends sensitively on the details of the physically underly- 
ing distribution functions, in particular the active black hole 
mass function and the Eddington ratio distribution. 

• Biases from measurement uncertainties: Both black hole 
mass and bulge property measurements suffer from measure- 
ment and calibration errors. In combination with selection 
limits they may introduce additional biases. 

• AGN evolution bias: The AGN distribution functions that 
regulate the magnitude of the bias are evolving themselves 
as a function of redshift, as implied by the evolution of the 
AGN luminosity function. This evolution will change the ex- 
pected sample bias with redshift. 

These effects alter the bivariate probability distribution and 
can bias the conclusions drawn from observational studies. It is 
however important to clarify what is meant by saying that mea- 
surements of the M.-bulge relations are biased. These relations 
can be represented in a number of ways that differ in complexity. 
(1) The most comprehensive approach would be to investigate 
the full bivariate distribution of galaxy properties and black hole 
masses. This distribution is modified by the AGN selection, and 
adopting the observed relation as the measurement thus consti- 
tutes a bias with respect to the true relation. (2) The most basic 
approach involves computing simply the mean offset of the en- 
tire sample from the local relation. We showed that this quantity 
is equally biased by the AGN selection. (3) An intermediate way 
is to look at the projected M.-bulge relation, defined as the mean 
black hole mass for a given galaxy property. This relation also 
suffers from selection biases, caused by excluding black holes 
with certain properties from the sample during the observational 
process. (4) Alternatively, one can evaluate the inverse projected 
relation, thus the mean galaxy property at a given black hole 
mass. While for our parameterisation of the bivariate distribu- 
tion this projection deviates from the true relation even without 
any selection effects, it is not modified by the AGN selection ef- 
fects outlined above. Therefore this relation has the potential to 
be used as an unbiased estimator for the M.-bulge relation, as 
long as additional selection on the host galaxy properties do not 
play a role. Among these four outlined approaches, options (2) 
and (3) are most commonly adopted in the literature. 

It would clearly be desirable if we could properly model and 
correct the observations for the known sample selection effects. 
In order to do so, it is necessary to know the relevant underly- 
ing distribution functions, i.e. the distribution of spheroid prop- 
erties such as velocity dispersions or masses; the active frac- 
tion as a function of BH mass, or alternatively the active black 
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hole mass function; and the Eddington ratio distribution func- 
tion. Thanks to recent observational progress, these distributions 
are now reasonably well known for the local universe; but they 
are still rather ill-constrained for higher redshifts. 

We employed our formalism together with a set of best-guess 
distribution functions to investigate how much existing studies of 
the M. -bulge relation in low-redshift AGN samples might suf- 
fer from selection biases. For the heterogeneous sample of AGN 
with reverberation mapping data we estimated a mean offset of 
~ 0.1 dex in M. at given bulge property, chiefly due to the nec- 
essary flux limitations imposed in the process of selecting these 
objects for observations. This is a small but not entirely negligi- 
ble effect, especially in view of the important role of this sample 
to serve as absolute calibration of the virial method for estimat- 
ing back hole masses. We also confirm that slope and scatter are 
not significantly affected. 

Considering the more general case of low-redshift AGN 
samples with measured bulge properties and BH masses esti- 
mated by the virial method, we find that the observed M, -bulge 
relations suffer from a strongly mass-dependent offset from the 
true relation. As a consequence, the observed relation of M. 
at given bulge property is flattened notica bly. Such a flattening 
has been noted several times by observers (iGreene & Hoi l2006t 



has been noted se veral times b y observers (ureene & Hot IzUUbt 
iKimetaU 120081: iBentz et al.L [2009bt IWoo et all 120101) : it can 



be explained in our framework without having to invoke intrin- 
sic differences in the BH - bulge relation between active and 
inactive galaxies. 

Extending this type of analysis to higher redshifts is ham- 
pered by our poor current knowledge of the various relevant 
distribution functions. In order to constrain the evolution in the 
M. -bulge relations it is crucial to not only obtain larger ob- 
served samples, but also to improve on measuring the under- 
lying distributions. For example, not knowing the mass depen- 
dence of the active fraction at high redshifts will fundamentally 
limit our ability to correct for the biases associated with this ef- 
fect. Furthermore, the selection function of the studied sample 
needs to be known; for an ill-defined heterogeneous sample the 
chances to properly account for the selection effects are heavily 
reduced. 

We briefly reviewed a number of recent observational at- 
tempts to study the M.-bulge evolution. The result is sobering: In 
no case do we find statistically significant evidence for an evolv- 
ing M.-bulge relation. While the observed apparent offsets in 
M. /Mfiuise from the local relation can be quite large, the sample 
selection bias estimated from our formalism is typically of the 
same magnitude. This does of course not exclude the possibil- 
ity of real evolutionary effects that might even be already visible 
in some of the data; but our ability to distinguish between real 
and artificial trends caused by selection biases is presently not 
sufficient. 

A possible route to circumvent several of the most problem- 
atic issues connected with AGN selection is the approach to fo- 
cus on the mean galaxy property at a given black hole mass. 
However, while this projected distribution is formally unaffected 
by the AGN-specific selection biases, it also does not directly 
yield the intrinsic relation. In order to reconstruct the true rela- 
tion in this case it is necessary to know the spheroid distribu- 
tion function and the intrinsic scatter in the M.-bulge relation. 
Neither is well constrained at high redshifts, but these quantities 
are observables, so that this approach may become very power- 
ful in the future. 
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Appendix A: Validation of the bivariate probability 
distribution 

Our chosen parameterisation of the bivariate probability distri- 
bution given as 



V(s,v)=g( M \s)<!> s (s), 



(A.l) 



(cf. Equation ([U), follows common practice, but it is not the 
only possible solution. In particular, this parameterisation pro- 
duces the visual impression of an upturn in the contours towards 
high values of s, corresponding to an upwards curved ridge line 
of the distribution. This upturn is simply a consequence of the 
fact that the space density of galaxies is steeply declining to- 
wards high masses or velocity dispersions. Since the quantita- 
tive conclusions in this paper depend on the appropriateness of 
the adopted form for T^s,//), it is worthwhile to test if our pa- 
rameterisation is actually consistent with current observations. 
Here we perform such a test. 

In Fig. lA.ll we superimpose our intrinsic bivariate probabil- 
ity distribution ^{s,^), which includes an estimate of the rele- 
vant bulge property distribution function, with data of galaxies 
for which dynamical black hole mass measurements are avail- 
able. We consider the distribution of M. versus the three most 
commonly measured galaxy properties: Stellar velocity disper- 
sions o%, K band absolute magnitudes M#, and stellar bulge 
masses Mb u i ge . For the M.-cr, relation we incorpora ted the stel- 
lar ve locity dispersion distribution function from ISheth et all 
and compare t he result to the observational sample from 
iGiiltekin et al.l (120091) . For the M.-M K relatio n we employed 
the 2M ASS K-band luminosity function f rom iKochanek et al.l 
(1200 ll) and show the sample from lHul (120091) for comparison, and 
for the M. -Mbuige relation we use our estimate for th e spheroid 



mass f unction, based on galaxy mass functions from Be ll et al 



( 2003 ). and compare this with the sample from lHaring & Rix 
(200#. 



It is important to note that in this comparison we intention- 
ally do not consider any sample selection effects. Thus, the dis- 
tribution of data points in the three diagrams does not, and can- 
not, agree with the overall probability distribution. The test that 
we want to perform focuses on the upper right corner of each 
diagram, and the question is: Do the data follow the mentioned 
upturn in the ridge line? While there is of course considerable 
observational and intrinsic scatter in the data, a visual inspection 
of Fig. lA.ll clearlv confirms that indeed there is such a trend. 

We can perform a more quantitative test of the adopted 
ansatz by separately looking at the two univariate projections 
of the bivariate distribution, {y){s) (Equation © - represented 
by the thick solid lines in Fig. IA.lt . and (s)(fx) (Equation (O - 
shown by the red dashed-dotted lines. To compare these curves 
with the measurements we computed the mean fi in bins of s and 
the mean s in bins of /i from the data, shown in Fig. IA. li as black 
open squares and red open triangles, respectively. Towards high 
values of s there is a clear trend for the triangles to be located 
above and to the left of the squares, and the predicted curves 
trace the points quite well. We performed a simple x 1 test on 
all points above a certain minimum value of s or fi, respectively, 
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Fig. A.l. Comparison of our bivariate probability distribution ^(s, fj) (contours) with observational data of black holes with dynam- 
ical mass measurements (blue circles). The thick black line shows the mean relation (fi)(s), the red dashed dotted line shows (s)(ji). 
The black open squares gives the mean fi binned in s from the observ ations (an estimate of (u)(s)) and the red open triangles gives 
s binned i n fj. (an estimate of (s}(ji)). Left panel: M. - tr, sample from l Giiltekin et all d2009l) . Middle panel: M. - Mk sample from 
iHul (120091) . Right panel: M. - M bu i ge sample from lHaring & Rixl (12004) . 



and found that the binned data are consistently described by the 
predicted curves, with a reduced x 1 always of order unity inde- 
pendently of the chosen cutoff. We conclude that our adopted 
form for the bivariate probability distribution ^(s,^) is consis- 
tent with the current observations. 



